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In this work we revisit the topic of two-dimensional Bose-Einstein condensates under the influence 
of time-dependent magnetic confinement and time-dependent scattering length. A moment approach 
reduces the examination of moments of the wavefunction (in particular, of its width) to an Ermakov- 
Pinney (EP) ordinary differential equation (ODE). We use the well-known structure of the solutions 
of this nonlinear ODE to "engineer" trapping and interatomic interaction conditions that lead to 
condensates dispersing, breathing or even collapsing. The advantage of the approach is that it is 
fully tractable analytically, in excellent agreement with our numerical observations. As an aside, we 
also discuss how similar time-dependent EP equations may arise in the description of anisotropic 
scalar field cosmologies. 

PACS numbers: 

I. INTRODUCTION 

The recent progress in experimental and theoretical studies of Bose-Einstein condensates (BECs) of dilute atomic 
gases [1] has been tremendous after their experimental realization Q. This has also led to an explosion of interest in 
the theme of nonlinear matter- waves such as dark [3[, bright m and gap [5[ solitons. Two-dimensional (2D) nonlinear 
excitations of BECs, such as vortices |6[ and vortex lattices [7], were also realized experimentally, while a multitude 
of other coherent nonlinear structures were also theoretically predicted; these include, e.g., Faraday waves [8], ring 
dark solitons and vortex necklaces [9], stable solitons and localized vortices in attractive BECs trapped in periodic 
potentials [10[, matter- wave gap vortices 2D bright solitons in dipolar BECs [12[, and so on. 

From the theoretical standpoint, the dynamics of such higher-dimensional structures is, generally, difficult to be 
treated analytically, especially so in time-dependent settings. This, in turn, raises important questions concerning their 
"controllability" , which would be of particular relevance regarding potential applications. Ideally, such a controllability 
would allow "manipulation" of the condensates at will, e.g., sustaining condensates which may expand towards 
complete dispersion, contract towards a complete collapse, or perform stable breathing oscillations. Moreover, such 
processes could involve a targeted growth of a condensate up to a certain width or shrinkage down a desired size. 
Our purpose in the present work is to illustrate how one can achieve this goal, by taking advantage of one of the few 
analytically tractable tools in higher-dimensional condensates, namely of the analysis of their moments [l3|, 03, Ql| 
(see also references therein). In particular, we focus on the case of quasi two-dimensional (2D), so-called, pancake- 
shaped condensates [l6| , under the effect of time- varying harmonic trapping pj], EH and also time-dependent s-wave 
scattering length (which controls the inter-particle interactions). Note that the controllability over the harmonic 
trapping is straightforwardly realizable under AC- variation of the atom trapping frequencies, while the controllabity 
over the interatomic interactions can be realized by using the so-called Feshbach resonance [17], connecting the s- 
wave scattering length to external magnetic fields. This latter technique (which is usually called Feshbach resonance 
management (FRM)) has been proposed theoretically as a means of avoiding collapse [18|, UJJ [20[ , but also as a way 
of producing robust coherent nonlinear matter- waves [2l|, [22|, [23|, |24| . 

For the quasi-2D condensates discussed above, upon presenting the moment analysis, we obtain a reduced dynamical 
description effectively involving only a variable associated with the width of the condensate wavefunction. This 
ordinary differential equation (ODE) can be reduced to one of the Ermakov-Pinney (EP) type [25|, [26|, [53, , whose 
solutions can be obtained analytically, provided that the solutions of the underlying linear Schrodinger equation can 
be obtained. We use this feature and the freedom in selecting the time-dependent trapping and interactions of the 
condensates to illustrate that one can construct analytical solutions to this ODE that completely characterize the 
temporal evolution of the width of the wavefunction. In so doing, we fully prescribe the dynamical evolution of 
the condensate. We show three prototypical examples in applying this idea: one in which the width grows in time 
(leading to condensate expansion), one in which it decreases (leading to focusing), and one in which it periodically 
"breathes" between a minimum and a maximum value. In all three cases, we test the analytical prediction against the 
full numerical simulation of the mean- field partial differential equation (PDE) model fully describing the condensate. 
We find excellent agreement between the two, showcasing the accuracy of our theoretical approach. 

The analysis of the ensuing EP equation in this setting is of interest in its own right. This is because, as we will 
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see, the resulting EP equation has a time- dependent nonlinearity in the right hand side (contrary to what is the case 
for the "standard" EP framework; see e.g. [25|, [26|, [2?], [28[ and references therein). In this context, and as an aside, we 
present a second physically relevant example where such EP equations with time-dependent nonlinearities may arise, 
by studying anisotropic scalar field cosmologies of a particular anisotropic geometry. This generalizes the example of 
(25^ parallelizing BECs without FRM and scalar field cosmologies in the isotropic case. 

Our presentation will be structured as follows: In section III, we provide a brief synopsis of the main features of the 
EP equation. In section IV, we proceed to an overview of the moment analysis, following the earlier works [l3L [13. fl^ . 
In section IV, we present the analytical solutions of the EP equation that we develop for each of the above-mentioned 
three cases. Section V tests these analytical results against full simulations of the PDE describing the condensate. 
Section VI illustrates the second example of time-dependent EP equations in the description of anisotropic scalar field 
cosmologies. Finally, in section VII, we summarize our findings and present our conclusions. 



II. THE ERMAKOV-PINNEY EQUATION 

The Ermakov-Pinney (EP) equation is a remarkable nonlinear ODE of the form: 

Y" + Q(r)Y=^. (1) 

The particularly attractive feature of this nonlinear ODE is that its general solution can be obtained, provided that one 
is able to solve the time-independent linear Schrodinger (LS) equation Y" + Q(r)Y = 0. For details on the properties 
of the EP equation, the interested reader is referred to [26, 27] and references therein. Here we just mention its basic 
superposition principle property. Namely, if the linearly independent solutions of the LS equation are Y\ (r) and Y 2 (r) , 
then the most general possible solution of the EP equation is given by 

Y{t) = (AY? + BYf + 2CY 1 Y 2 ) 1/2 (2) 
where A, B and C are constants connected through 

AB ~ c2 = w (3) 

where W = Y{YJ> - Y 2 Y{ is the Wronskian of Yi(r) and Y 2 {r). 



III. MOMENT ANALYSIS FOR BECS 



One of the interesting variants of the "regular" EP equation of the form ([T]) arises in the study of BECs, albeit in 
a somewhat modified form (see below). 

The relevant mean-field model for studying atomic Bose-Einstein condensates at zero temperature consists of the 
so-called Gross-Pitaevskii equation [1] of the following dimensionless form: 

id t u = -^Au + (A(t)r 2 + v(t)\u\ 2 ) u, (4) 

where u represents the wavefunction of the condensate, V(r) = X(t)r 2 denotes the harmonic trap confining the bosons 
[with a time-dependent frequency determined by X(t)] and v(t) is the coefficient of the nonlinear term proportional 
to the s-wave scattering length, characterizing the interaction between the particles. We will take advantage of the 
Feshbach resonance [17j to consider that the latter is time-dependent, as well. 

One of the popular approaches to studying Eq. (j4]) is through the use of moment methods [H, Q3, [H| . The latter 
allow us to write ODEs for the moments of the mean-field wavefunction u as follows. We define 
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where subscripts denote partial differentiation, denotes complex conjugate and (d) indexes the dimension (similarly 
to IHIiil). The Hamiltonian of Eq. Q 

H=^ [\X7u\ 2 + v(t)\u\ A + 2A(t)r 2 |^| 2 ] r d ~ l dr, (9) 

2 .In 



can then be written as: 
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One can then infer from the dynamics of Eq. (J2J) that: 



\v{t)i$ + \{t)M. (io) 



tf=^$ + A'(*)J$ (11) 

However, as derived in pH, HEj. 

/g = 4Jf - 8A(t)I$ + 2i/(t)(d - 2)4$. (12) 

We focus on , since this moment is associated with the width of the spatial profile of the wavefunction. From the 
latter, we can infer a number of useful pieces of information concerning certain asymptotic values; in particular, if 
1^2 ~~ ^ 0? tlie- condensate collapses, while if —> oo, the BEC disperses. Furthermore, since the L 2 norm of the 
wavefunction is conserved, the estimate of on the (square) width of the wavefunction can be used together with 
this conservation law to provide information on the amplitude A of the wavefunction (i.e., approximately A ~ 1//^). 
Combining the two equations (jTTJ) and ([12]) . we obtain a single equation for the time-dependence of as 

lf 2 = 2i/(t)4J - 4A'(*)I$ - 8A(t) /g + (d - 2) [2i/(t)I$ + M*) (13) 

Following [gjj, we consider a quadratic phase for the solution (an assumption most relevant to dimension d — 2 as 
discussed in [20[, but which we will also consider more generally), and obtain (cf. Eq. (6d) of |20|) that 

I {d) - — (U) 

i 5,o - ( d ) ' v iZt ; 

1 2,2 

where the constant K is determined by initial conditions. Notice that this is the only assumption in our calculations 
herein, whose validity will be examined a posteriori by comparing our analytical results with numerical computations. 
Denoting for simplicity = y, the resulting ordinary differential equation for y can be written as: 

| (yy - i*> - + W),>) = m* - %| (f) <u> 

Clearly, from the above exposition, the most straightforward case is the one with d = 2, on which we will focus 
next. When d = 2, the equation can be directly integrated, yielding 

.. 1 . 2 . , > 2Kv(t)+C 

V ~ Yy V + mV = 7 ' ( 6) 

where C is an integration constant that can be computed from the initial conditions as: 



5(0) (if (0) - A(0)lg(0)) - l(lgl(0)) 2 - 2Kv(0) (17) 



C = 4lP 

If we now use the transformation Y(t) = y/y(t) \l4. l25j|. then we obtain an Ermakov-Pinney (EP) type equation of 
the form: 

Y + 2\{t)Y= Kv{t ^^ (18) 
Y 6 

Hence, Eq. ([T8]) is the equation that describes the dynamics of two-dimensional Bose-Einstein condensates in the 
presence of a time-dependent trap [H, [3, HBf , as well as in the scenario of Feshbach resonance management [H, fl9l . 

[Si, mi, m, . 

One can now examine particular cases of time-dependence of \{t) and ^(t), using the analytical tract ability of the 

(2) 

ensuing EP equation, in order to obtain completely analytical solutions for Y(t), and hence for I\ 2 . 
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IV. ANALYTICAL RESULTS 

There are numerous possiblities in the case of a time-dependent A and thus, as explained above, we limited our 
investigation to the three cases of an expanding waveform, a collapsing waveform and an oscillatory waveform. In 
order to derive examples for each of these cases, v(t) was chosen to be independent of the choice of X(t) and y(t). 
Solutions to Eqn. ([T8]) were determined by "reverse engineering" [28[ (based on the desired behavior of the dynamics 
of the wavefunction) , and the details of the EP functional form. Below are the three cases examined. Notice that in 
all three cases, the function f(t) used below is given by f(t) = Kv(t) + C/2 i.e., by the numerator of the right hand 
side of Eq. (fTHj) . 

1. The first case explores the possibility of an expanding width wavefunction. In order for the wavefunction width 
to increase, X(t) must be decreasing, thus test functions for X(t) and y(t) were chosen to reflect this. A similar 
methodology was also used in the two other cases to determine test functions. 



y(t) = B 2 (A 2 +t 2 ) (19) 
= f(t)-A*B* 

V 2£ 4 (A 2 + t 2 f V ) 

z(t) = B A A 2 = const. (21) 



2. The second case, a focusing wavefunction, required an additional condition (y(0) ^ oo) in order to perform 
the numerical simulations described below. Additionally, in this case and the next, an initially undetermined 
function z(t) was used in X(t) to simplify the derivation of X(t) and y(t). After these were found, the exact form 
of z(t) was calculated by plugging \(t) and y(t) into Eqn. ([T8|) and solving for z(t). This results in 

y(t) = ( 22 ) 



A 2 +tP 

X(t) = [f(t)-z(t)] B " (23) 
[2(A 2 +V>)] ql 



z(t) = 
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3. Finally, in the last case, equations for X(t) and y(t) were calculated in order to generate a wavefunction with 
oscillatory width. The corresponding functions in this case read 

y(t) = 6 + sin 2 (ct) (25) 

m . r f» -'W (26, 

2 [6 + sin 2 (ct)] 
c 2 

z(t) = -— sin 2 (2ct) +c 2 cos(2ct) [6 + sin 2 (ct)j (27) 

The resulting form of the normalized confining frequency X(t) for the three different cases is shown in Fig. [U 
We note in passing that such schemes are definitely realizable within the unprecedented control that exists over the 
magnetic confinement of the condensates. A very recent example illustrating this point can be found in the recent 
work of [29[; this experiment clearly realizes a form of the third type of confinement (an oscillatory one), by using 
a periodic modulation of the transverse confinement of the condensate, in order to induce longitudinal oscillations 
and the emergence of Faraday patterns. This directly shows that temporal modulation of the trapping frequencies is 
feasible. What we propose here consists of the three principal types of trapping frequency dependence, i.e., monotonic 
decrease, monotonic increase and oscillatory dependence, that allow an explicit analytical handle on the dynamics of 
the BEC, in excellent agreement with our numerical findings (see below). Hence, we expect that such schemes would 
be directly applicable in experimental settings similar to those of the above experiment. 



V. NUMERICAL SETUP AND RESULTS 



In order to validate our analytical results based on the predictions made above, the evolution of the condensates 
under the time-dependent trappings imposed by Eqs. ([20]) , ([23]) and ([26]) were tested by numerically solving the 



Time(t) 



Time(t) 



Time(t) 



FIG. 1: Plots of X(t) for Cases 1-3, i.e., respectively expanding (left), contracting (middle) and oscillatory (right) condensates. 

original PDE model of Eq. The numerical simulations were run using spectral methods [30[, with a spatial 

grid composed of Chebychev nodes (x n = cos(7rn/N) where N is the total number of nodes). Since the Chebychev 
nodes are contained within the interval [—1,1], the spatial variable was normalized in order for the solution u(t) to 
be supported on the interval [—1,1] from t = to some evolution horizon t = T > 0. The temporal integration was 
implemented using a 4th-order Runga-Kutta scheme. 

It is worth noting here that in all three cases examined, the choice of v(t) was independent of the equations for 
X(t) and y(t), thus v(t) = sin 2 (£) was chosen to be used for all cases, in consonance with the Feshbach resonance 
management scheme discussed in the Introduction; in that sense, however, notice that the key temporal variation in 
the schemes presented herein is that of the magnetic trapping frequency. As an initial condition to Eqn. (|4|), we used 
a generic Gaussian profile of the form u(t) = Aexp(— Br 2 /2) with the values of A and B subject to the condition 
y(0) = A 2 /(2B 2 ). 

For the cases of the expanding and collapsing width wavefunctions, the numerical results are in perfect agreement 
with the analytical calculations as shown in Figs. ([2|) and ([3]). Notice in the lower left-hand plot for both figures how 
the lines for the analytical and numerical plots are completely indistinguishable from each other. 

For the last case of the oscillatory width wavefunction, there was a very small discrepancy between the analytical 
calculations and numerical simulation, as seen in the lower left-hand plot of Fig. (|4|). This discrepancy is only seen 
at the maximum values of the wavefunction width, starting with the second and subsequent maximums. While this 
small discrepancy may be triggered by the sole approximation of our approach, namely Eq. (Q3j), further numerical 
experiments (not shown here) seem to indicate that it is more likely to be the result of the numerical approximation 
to the corresponding moment. 

(2) 

In any case, the overall excellent numerical agreement between the analytically obtained moment I\ 2 anc ^ ^ s 
numerically found counterpart illustrate the relevance of this approach and its usefulness in systematically prescribing 
the wavefunction behavior, in an analytically tractable way. 

VI. ANISOTROPIC SCALAR FIELD COSMOLOGIES 

An interesting analogy has been recently made between Bose-Einstein condensates and isotropic scalar field cos- 
mologies in the absence of time dependence in the coefficient of the nonlinearit y 12511. The so-called Friedmann- 
Robertson- Walker (FRW) metric, coupled to a scalar field had been shown earlier [26, 28] to be described by an EP 
equation and that viewpoint was exploited in [25[ to illustrate an analogy between cosmic dynamics and BECs. 

Here, we extend this analogy to EP equations with time-dependent nonlinear it ies. In particular, we show that 
anisotropic cosmic dynamics can be parallelized with Feshbach-managed Bose condensates. 

If we consider the spatially homogeneous, yet anisotropic geometry of Bianchi type I, we have a line element of: 

ds 2 = -N(t) 2 dt 2 + A(t) 2 dx 2 + B(t) 2 dy 2 + T(t) 2 dz 2 . (28) 

We can define the tensor: F^ v — G^ v — SnT^ where G^ v = R^ v — ^g^ u R is the Einstein tensor and T^ v = 
0;/x</V — \9^v{ < i>' (X( i>\a, -\-m 2 ()) 2 ) is the energy momentum tensor. Then, the quadratic constraint is the equation F$ = 0, 
the kinematic equation is given by Fl =0, while the Klein-Gordon equation for the scalar field (that is coupled to 
gravity) is given by — m 2 (j) = ex T.^ v = 0. Notice additionally, that the two integrals of the motion, namely 
h = Fl - F 2 = and J 2 = Fl - F 3 3 = yield B(t) = A(t)e Kt l 2 and T(t) = A(t)e xt / 2 . 
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FIG. 2: Evolution plots of the expanding width wavefunction of Eqs. (fT9]) - ([20|) . The upper left figure shows the spatio- 
temporal evolution of a fixed contour of the solution over time. The lower left figure compares the analytically determined 
moment associated with the wavefunction width with its numerically computed counterpart. The remaining figures show the 
contour of the solution at various times during the evolution, clearly indicating the expanding nature of the BEC wavefunction. 



Solving the Klein-Gordon equation for <j)"(t) and substituting the result into dtF§ = (as well as solving F® = for 
(j) and using the resulting expression in d t F® = 0), one is led to a dynamical equation for the remaining scale factor 
A(t) in the form: 



KA'jt) \A\t) A A' (if ftff A"(t) 
4 + A(t) + A(t) + A {tf 2 A(t) 

Using now: A(t) = Y(t) 2 / n and a change of variable r = J Sl(t')dt' , we obtain: 



Y(t) + Y(t) 



Q'(t) r(r) 2 (6 + n) 



'm 2 y{t) 

Hence, a choice of time reparametrization according to: 



(k + A) 



Y(t) uk\ Y(t) nY{ T )4>{T) 2 



n(t) 2 



(6 + n) Y'(t) 



Q(t) n Y{t) 

(which leads to fl(t) = 6»e( K+A ^y(t)( 6+ra )/ n , where 9 > is a constant of integration), results in the form: 



Y(t)+Q(t)Y(t) 



T(t(r)) 

y( T )l + 12/r, 



>(r) 2 



8^2 ' 



(29) 



(30) 



(31) 



(32) 



Eq. ([32]) becomes an EP equation for the choice of n = 6. Notice that in Eq. ([32]) , a nontrivial complication is that the 
time dependent coefficient T depends on the reparametrization of time through Q. Nevertheless, the original equation 
([29]) has been solved in [3l|. This may, in turn, provide valuable insights in the solution of such time-dependent EP 
equations; this would be an interesting direction for future studies in such time-dependent EP settings. 

We close this section by noting that one can therefore generalize the analogy of Bose condensates and scalar field 
cosmologies in the Feshbach-managed viz. anisotropic case of the EP equation with time-dependent nonlinearity. The 
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FIG. 3: Evolution plots of the focusing case of Eqs. ([22 ]) -([24 ]) . with p = 2 and q = —9. The upper left figure shows the spatio- 
temporal evolution of a fixed contour of the solution over time. The lower left figure compares the analytically determined 
wavefunction width with its numerically computed analog. The remaining figures show the contour of the solution at various 
times, illustrating the focusing effect. 

quantities that now bear the analogy is the scale factor of Bianchi type I with the second moment of the condensate 
wavefunction, the time-dependent magnetic trap strength with the time-dependent scalar field and finally the time 
dependent nonlinearity coefficient (i.e., the scattering length in the condensate dynamics) is (indirectly) connected 
with the reparametrization of time in the cosmological problem. 

VII. CONCLUSIONS AND FUTURE CHALLENGES 

In this short communication, we revisited the theme of higher dimensional Bose-Einstein condensates under the 
presence of magnetic trapping and Feshbach Resonance management. We used the moment method to develop an 
Ermakov-Pinney ODE for the second moment of the distribution function, which is associated with the width of the 
condensate. This EP equation is analytically tractable, in a number of cases. In particular, one can "reverse engineer" 
magnetic trappings that will induce the expansion, contraction or oscillatory behavior of the condensate at will. For 
such scenaria, the relevant moment of the wavefunction can be obtained analytically and is found to be in excellent 
agreement with our numerical simulations of the full original PDE model. The approach permits a detailed analytical 
handle on the behavior of higher dimensional condensates which is usually quite difficult to acquire with different 
methods based on nonlinear PDEs. It may also, in turn, permit to appropriately craft experiments, based on the 
reshaping "operation" that is desirable to perform on the condensate. 

As an aside example, we have illustrated the relevance of such time-dependent Ermakov-Pinney equations in a 
completely different physical system, namely in anisotropic scalar field cosmologies of Bianchi type I. 

This approach also suggests a number of interesting questions. It would be, in particular, relevant to examine 
whether a general solution of the time-dependent EP equation developed herein can be obtained on the basis of 
its linear Schrodinger counterpart. This would bear direct consequences both in the atomic physics and in the 
cosmological problem, allowing for a general analytical description of their respective time-dependent properties (for 
the condensate wavefunction width or the cosmological scale factor, respectively). Furthermore, the approach was 
presented here for two-dimensional settings for reasons that have to do with the simplification/closure of the moment 
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FIG. 4: Evolution plots for the oscillatory width wavefunction case of Eqs. (|25]) - (|27j) . The panels are similar to the previous two 
figures. The middle and right panels show the contour of the solution at the first two maxima and minima of the wavefunction 
width. 

approach in that case. However, it would be of particular interest to develop similar approaches and potential closure 
schemes for one-dimensional or three-dimensional settings. Such studies are currently in progress and will be reported 
in future publications. 



[1] F. Dalfovo et a/., Rev. Mod. Phys. 71, 463 (1999). 

[2] M.H. Anderson et al, Science 269, 198 (1995); K.B. Davis et al, Phys. Rev. Lett. 75, 3969 (1995); C.C. Bradley et al, 
Phys. Rev. Lett. 75, 1687 (1995); see also the reviews E. A. Cornell and C. E. Weiman, Rev. Mod. Phys. 74, 875 (2002) 
and W. Ketterle, Rev. Mod. Phys. 74, 1131 (2002). 

[3] S. Burger et al, Phys. Rev. Lett. 83, 5198 (1999); J. Denschlag et al, Science 287, 97 (2000); B.P. Anderson et al, Phys. 
Rev. Lett. 86, 2926 (2001). 

[4] K.E. Strecker et al, Nature 417, 150 (2002); L. Khaykovich et al, Science 296, 1290 (2002). 

[5] B. Eiermann et al, Phys. Rev. Lett. 92, 230401 (2004). 

[6] M.R. Matthews et al, Phys. Rev. Lett. 83, 2498 (1999); K.W. Madison, et al Phys. Rev. Lett. 84, 806 (2000); S. Inouye 

et al, Phys. Rev. Lett. 87, 080402 (2001). 
[7] J.R. Abo-Shaeer et al, Science 292, 476 (2001); J.R. Abo-Shaeer, C. Raman and W. Ketterle, Phys. Rev. Lett. 88, 070409 

(2002); P. Engels et al, Phys. Rev. Lett. 89, 100403 (2002). 
[8] K. Staliunas, S. Longhi and G.J. de Valcarcel, Phys. Rev. Lett. 89, 210406 (2002). 
[9] G. Theocharis et al, Phys. Rev. Lett. 90, 120403 (2003). 
[10] B.B. Baizakov, B.A. Malomed, and M. Salerno, Europhys. Lett. 63, 642 (2003); D. Mihalache et al, Phys. Rev. A 72, 

021601 (2005). 

[11] E.A. Ostrovskaya and Yu.S. Kivshar, Phys. Rev. Lett. 93, 160405 (2004). 
[12] P. Pedri and L. Santos, Phys. Rev. Lett. 95, 200404 (2005). 

[13] J.J. Garcia-Ripoll and V.M. Perez- Garcia, [patFsol/9904006 ; see also V.M. Perez-Garcfa et al, J. Opt. B 2, 353 (2000). 

[14] J.J. Garcia-Ripoll, V.M. Perez-Garcfa, and P. Torres, Phys. Rev. Lett. 83, 1715 (1999). 

[15] P.G. Kevrekidis, A.R. Bishop and K.0. Rasmussen, J. Low Temp. Phys. 120, 205-212 (2000). 

[16] A. Gorlitz et al, Phys. Rev. Lett. 87, 130402 (2001); M. Greiner et al, Phys. Rev. Lett. 87, 160405 (2001); S. Burger et 
al, Europhys. Lett. 57, 1 (2002). 



9 



[17] S. Inouye et a/., Nature 392, 151 (1998); J. L. Roberts et a/., Phys. Rev. Lett. 81, 5109 (1998); J. Stenger et a/., Phys. 

Rev. Lett. 82, 2422 (1999); S. L. Cornish et a/., Phys. Rev. Lett. 85, 1795 (2000). 

[18] F.Kh. Abdullaev et a/., Phys. Rev. A 67, 013605 (2003); H. Saito and M. Ueda, Phys. Rev. Lett. 90, 040403 (2003). 

[19] F. Kh. Abdullaev et a/., Phys. Rev. A 68, 053606 (2003). 

[20] G.D. Montesinos, V.M. Perez-Garcfa, and P.J. Torres, Physica D 191, 93 (2004). 

[21] P.G. Kevrekidis et a/., Phys. Rev. Lett. 90, 230401 (2003). 

[22] D.E. Pelinovsky, P.G. Kevrekidis and D.J. Frantzeskakis, Phys. Rev. Lett. 91, 240201 (2003). 

[23] Z. Rapti et a/., Phys. Scripta T107, 27-31 (2004). 

[24] F.Kh. Abdullaev et al, Phys. Rev. Lett. 90, 230402 (2003). 

[25] J.E. Lidsey, Class. Quantum Grav. 21, 777 (2004). 

[26] R.M. Hawkins and J.E. Lid sey, Phys. Rev . D 66, 023523 (2002). 

[27] P.B. Espinoza Padilla, [math-ph/ 0002005 

[28] F.L. Williams and P.G. Kevrekidis, Class Quantum Grav. 2 0, L177 (2003). 

[29] P. Engels, C. Atherton and M.A. Hoefer, [cond-mat/07010281 

[30] L.N. Trefethen, Spectral Methods in MATLAB, SIAM (Philadelphia, 2000). 

[31] T. Christodoulakis, Th. Grammenos, Ch. Helias, P.G. Kevrekidis and A. Spanou, J. Math. Phys. 47, 042505 (2006). 



